/**
 * TODO: Should probably NOT modify the Ising input class. This can be done
 * by creating a new Ising instance, or perhaps by just being an extension of
 * Ising (inheritance).
 */

import java.lang.System;

public class Generate // Code to generate droplets
{
	final int numTests  =  100;
	final int upRelax   =   10;

	// nucCheck must be way less than Ising nucleation time
	final int nucCheck  =  150; // Number of steps to allow m to become < thres
	final double thresh = 0.75; // The nuc. has already occurred threshold

	// The amount of time possibilities the log-generate function searches
	final int dropTimeSearch = 128;

	Ising ic;
	long theSeed;

	public Generate (Ising ic) { this (ic, System.currentTimeMillis ()); }

	public Generate (Ising ic, long theSeed)
	{
		this.ic = ic;
		this.theSeed = theSeed;
	}

	// A very rough initial guess: nucleation when magnetization crosses 0,5
	private int initialGuess ()
	{
		int t = 0;
		ic.setSeed (theSeed);
		ic.allSpinsUp (); ic.upH ();
		for (int i=0; i<upRelax; ++i) ic.update ();
		ic.flipH ();
		while (ic.magnetization () > thresh) { ++t; ic.update (); }
		return t;
	}

	/**
	 * Recursively configure the internal Ising class to have a nuc droplet.
	 * Start with t above t_n, and once t below t_n with high enough confidence
	 * then we call this time t := t_n. The hypothesis test is for p greater 
	 * than 0.5.
	 * @param t the last nucleation time guess
	 */
	private void genDroplet (int t)
	{
		int numNucleate = 0;

		// Prepare the ising class
		ic.setSeed (theSeed);
		ic.allSpinsUp (); ic.upH ();
		for (int i=0; i<upRelax; ++i) ic.update ();
		ic.flipH ();

		for (int i=0; i<t; ++i) ic.update (); // Regenerate the state
		boolean [] s = ic.getLattice ();

		ic.randomizeSeed ();
		for (int i=0; i<numTests; ++i)
		{
			for (int j=0; j<ic.N; ++j) { ic.setSpin (j, s[j] ? 1 : -1); }
			for (int j=0; j<nucCheck; ++j) 
			{
				ic.update ();
				if (ic.magnetization () < thresh) break;
			}
			if (ic.magnetization () < thresh) ++numNucleate;
		}
		double pbar = (double) numNucleate/numTests;
		double sigma = 0.5; // Standard deviation for p = 1/2
		double z_alpha = 1.65; // 95% confidence (one sided)
		
		System.out.println (t+" "+pbar);
		if (pbar > 0.5 + z_alpha*sigma/Math.sqrt (numTests)) 
			genDroplet (t-1);
		else // restore droplet
			for (int j=0; j<ic.N; ++j) { ic.setSpin (j, s[j] ? 1 : -1); } 
	}
	
	/**
	 * Recursively configure the internal Ising class to have a nuc droplet.
	 * Same as other genDroplet, but performs in logarithmic time
	 * @param t the last nucleation time guess
	 */
	private void genDropletLog (int t, int search, int twoPow)
	{
		int numNucleate = 0;
		final double limit = 0.1;

		// Prepare the ising class
		ic.setSeed (theSeed);
		ic.allSpinsUp (); ic.upH ();
		for (int i=0; i<upRelax; ++i) ic.update ();
		ic.flipH ();

		for (int i=0; i<t; ++i) ic.update (); // Regenerate the state
		boolean [] s = ic.getLattice ();

		ic.randomizeSeed ();
		for (int i=0; i<numTests; ++i)
		{
			for (int j=0; j<ic.N; ++j) { ic.setSpin (j, s[j] ? 1 : -1); }
			for (int j=0; j<nucCheck; ++j) 
			{
				ic.update ();
				if (ic.magnetization () < thresh) break;
			}
			if (ic.magnetization () < thresh) ++numNucleate;
		}
		double pbar = (double) numNucleate/numTests; // estimator
		double sigma = 0.5; // Standard deviation for p = 1/2
		double z_alpha = 1.96; // 95% confidence (two sided)
		
		//System.out.println ("t="+t+", p="+pbar);
		
		// in case we can't land close enough to 0.5 (keep going above and 
		// below the nucleation limit), also only accept mature droplets
		if ((double) search/twoPow < limit && pbar > 0.5) // restore droplet
		{
			for (int j=0; j<ic.N; ++j) { ic.setSpin (j, s[j] ? 1 : -1); } 
			return;
		}

		if (pbar > 0.5 + z_alpha*sigma/Math.sqrt (numTests)) 
			// Extra 1 is to make sure we don't get stuck subtracting zero
			genDropletLog (t-search/twoPow-1, search, twoPow*2);
		else if (pbar < 0.5 - z_alpha*sigma/Math.sqrt (numTests)) 
			genDropletLog (t+search/twoPow+1, search, twoPow*2);
		else // restore droplet
			for (int j=0; j<ic.N; ++j) { ic.setSpin (j, s[j] ? 1 : -1); } 
	}

	public void genDroplet () 
	{ 
		//genDroplet (searchGuess ()); 
		int t = initialGuess ();
		genDropletLog (t, dropTimeSearch, 2);
	}
}


